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We show that there exists a long-range RVB state for the kagome lattice spin-1/2 Heisenberg 
antiferromagnet for which the spinons have a massless Dirac spectrum. By considering various 
perturbations of the RVB state which give mass to the fermions by breaking a symmetry, we are 
able to describe a wide-ranging class of known states on the kagome lattice, including spin-Peierls 
solid and chiral spin liquid states. Using an RG treatment of fluctuations about the RVB state, we 
propose yet a different symmetry breaking pattern and show how collective excitations about this 
state account for the gapless singlet modes seen experimentally and numerically. We make further 
comparison with numerics for Chern numbers, dimer-dimer correlation functions, the triplet gap, 
and other quantities. To accomplish these calculations, we propose a variant of the SU(N) theory 
which enables us to include many of the effects of Gutzwiller projection at the mean-field level. 
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I. INTRODUCTION 



The spin 1/2 Heisenberg antiferromagnet on the kagome lattice is a good candidate for a two-dimensional quantum 
system with a spin disordered ground state Q. While it appears that on square j| and triangular |||-||] lattices an 
antifcrrmomagnet will acquire Neel order, on the kagome lattice strong numerical evidence has accumulated that 
the system is spin disordered, as seen by the existence of a gap to triplet excitations and through consideration of 
the spectra of finite size samples || . Numerically, one finds a continuum of low energy states below the triplet gap 
. The continuum of low energy excitations provides a great puzzle to theory in the absence of an obvious broken 
O \ symmetry. 

i^i - There are good experimental realizations of kagome systems, despite the presence of additional couplings, including 
the jarosites and SrCrGaO. While in iron jarosites H|), these additional couplings produce long-range order, in 
deuteronium jarosite and SrCrGaO no long range order is seen. Additionally, in SrCrGaO a quadratic 
specific heat and very weak field dependence of the specific heat are in agreement with the picture of a continuum 
of low-energy singlet excitations seen in numerics, suggesting that the latter two compounds provide good realizations 
' of the kagome antiferromagnet. 
ly-j \ Given the lack of spin order, RVB ideas seem natural for this system, and indeed have stimulated much theoretical 
work on the system. Large N calculations based on SU(N) have been used to suggest a spin-Peierls state [ p^[ . 
Calculations based on Sp(2N) have suggested a phase with deconfined, gapped, bosonic spinons |l4|]. Chiral states 
have also been proposed [ fl3[ , but do not account for the excitation spectrum and also are in disagreement with the 
rapid decay of chirality-chirality correlation functions seen in numerics. States with BCS pairing have been suggested 
but again do not account for the excitation spectrum; due to the non-bipartite nature of the kagome lattice, these 
states are not equivalent to flux states [fL5[ - In addition to the long-range states, short-range RVB states based on 
a reduced Hilbert space of dimers |lq-|l8| have also been considered and provide some explanation for the gapless 



continuum. 



An RVB state on the kagome lattice would be particularly attractive, given the intensive work on RVB states on the 
square lattice [|l9| , especially in connection with high-T c materials J2(| ■ In the absence of doping, the square magnet 
eventually acquires Neel order and the spinon excitations disappear from the system. Since the kagome lattice does 
not acquire Neel order, it could be a very important model for a spin liquid or spin solid state. 

The idea behind the present approach is a to start with a long-range RVB treatment of the kagome lattice and 
consider various ways of gapping the spinon excitation spectrum. We will first construct a "parent state" which will 
be the best RVB state that does not break time-reversal symmetry or any lattice symmetry. 

We will then demonstrate an interesting massless Dirac structure for this state. Various other known RVB states 
can be obtained by perturbing the parent state, lowering the symmetry and giving mass to the Dirac particles, so 
that the parent state unifies a wide class of states. Physically, we expect that the system will attempt to give mass 
to the Dirac particles and open a gap, picking out one of these other lower symmetry states. We will discuss the 
symmetry breaking through a renormalization group treatment. We will obtain some kind of spin solid state, and 
some low energy Goldstone and gauge excitations, which we will argue provide the low energy degrees of freedom seen 
experimentally. 
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The RVB states can be thought of by a decoupling procedure, in which we decompose spin-1 operators into pairs 
of spin-1/2 operators. Take a Hamiltonian 

H = ^2 J Si ■ Sj (1) 

<i,3> 

where the sum extends over neighboring sites 

Introduce the spinon fields ip^i), ipa(i), where a — u,d labels up and down spinon fields. We can then introduce a 
Hubbard-Stratonovich held 4y = such that 

H = i^l^Mi) + h.c) ( 2 ) 

<i,3> <i,j> 

By taking a mean-field in t, minimizing the total energy of the fermions and the Hubbard-Stratonovich field, we obtain 
an RVB state. One must at some point project results onto the physical space in which each site is singly occupied. 

Later, we will find this projection to be extremely important. In the absence of projection, the ideal mean-field 
state is almost always found by taking a dimer covering of the lattice pi] , with nonvanishing only on the given 
dimers. Projection can stabilize RVB states, so although our first calculations will ignore the effects of projection, in 
a naive mean-field, we will later discuss a projected mean-field that includes some of the essential effects of projection. 

We will then have to proceed beyond mean-field solutions. We will consider a functional integral with fields tp(i) and 
Uj , fluctuating about a saddle-point of the action. There are a large number of possible fluctuations in tij , including a 
set of pure gauge fluctuations, as well as a set of gauge fields. Most other fluctuations can be ignored because they do 
not contribute to the low-energy dynamics. However, there will be a particular set of fluctuations in ty that produce 
a mass for the fermion field. Although these fluctuations are not gapless, we will retain these fluctuations due to their 
impact on the low energy dynamics of the fermion field. We will see using a renormalization group that the effective 
action of these fields can differ greatly from that suggested by the mean-field. 

To outline the paper, we will first describe the parent state, and then discuss how to perturb the parent state to 
obtain other proposed RVB states. Then we will discuss naive and projected mean-field theory treatments of these 
states. We will the proceed to a field theoretic treatment of fluctuations about the mean-field and a renormalization 
group that will suggest one particular symmetry breaking pattern. We will discuss the pseudo-Goldstone and gauge 
modes that arise from this symmetry breaking, and the mechanism that ultimately gives them a very small energy 
gap. 

Next, we proceed to a discussion of finite system size effects as a first step in comparison with numerics. These 
effects lead to an additional flux for odd system sizes which leads to a nonvanishing Chern number for odd system 
sizes. We will then compare the low energy bosonic modes from the field theory to the low energy singlet modes found 
in numerical calculations as well as checking dimer-dimer correlation functions and many-body density of states. 



II. THE PARENT RVB STATE 



Although the short range RVB calculations provide one starting point for the kagome lattice antiferromagnet, we 
will be interested in looking at long range states instead. Certainly, the short range RVB calculations themselves 
suggest that long range antiferromagnetic correlations are important; the variational energy of these states improve 
when second neighbor dimers are included. Further, while the kagome lattice has a gap to triplet excitations, this 
gap is about an order of magnitude smaller than J; from a short range calculation one might expect a gap of order J 
as that is the energy to break a dimer. However, Mila has suggested that within a short-range state it is possible to 
have a triplet gap much smaller than J jig]] , so the small triplet gap does not necessarily rule against a short-range 
state. 

The best RVB state on the kagome lattice antiferromagnet is a chiral spin liquid j|3]. A similar chiral state |p2 23 



was obtained using a hardcore boson representation of the spins and transmuting the statistics from bosonic to 
fermionic using a Chern-Simons field. However, numerical calculations p4| do not support a large chirality-chirality 
correlation function or expectation value of the chirality operator, which would seem to rule these states out. So, we 
will look for the best RVB state that is not chiral. 

Assuming that we are looking for a long-range RVB state in which all £y have the same magnitude, the only choice 
we have is how much flux to put into the system. The state we choose involves putting n flux through the hexagons, 
and no flux through the triangles. This state offers a better mean field energy than any other nonchiral RVB state, 
including the state with no flux through the system at all. 
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The unit cell of the kagome lattice consists of three sites on a triangle. Once we add flux to the system, the unit 
cell doubles, and requires six sites on two triangles. We will find it convenient to double the unit cell again, to twelve 
sites, including six sites on a hexagon and the six sites which neighbor the hexagon. This cell is shown in figure 1. 
The kagome lattice is made up of a triangular lattice of these 12-site unit cells. We have numbered the points in the 
cell for later reference. 

For now, let us assume that we pick J such that \t\ = 1 within our RVB state. Then, the band structure for our 
RVB state is shown in figure 2 scanning along the given line of momenta in the Brillouin zone. There is a degeneracy 
of states: the bottom line in the figure actually consists of four bands, while the other four lines in the figure consist 
of two bands each, providing a total of twelve bands. At (0, 0), four bands meet at energy less than zero and another 
four meet at positive energy. Near this point the spectrum becomes relativistic. 

The particles occupy the lowest six bands of the system, meaning that where the bands meet the spectrum becomes 
gapless. The system can gain energy by perturbing about our given RVB state; we expect that the greatest gain in 
energy comes from opening a gap. For this reason, we will study the Dirac point and look at possible perturbations 
to the Dirac equation. 

At the Dirac point, the Schrodinger equation for the fermions becomes 

Eip — Vf(a x k x + a y k y )ip (3) 

where ip is a four component spinor, and the matrices a x ,a y are anti-commuting a-matrices. The particular basis 
chosen for ip and for the a- matrices is unimportant. We find by explicit computation that 

v f = (0.408248.. .)\t\ (4) 

Given the Dirac equation (^), we would like to consider the effect of perturbations Stij on the low energy structure. 
This analysis will enable us to focus on those fluctuations in tij which have the greatest impact on the low-energy 
dynamics and which must be kept when we proceed to a field-theory treatment of fluctuations. 

From equation (j|), the system pays an energy cost equal to j\5tij\ 2 , but it can gain energy by opening a gap for 
the Dirac particles. As a result, we look for the perturbations which open the greatest gap for the Dirac particles for 
given |(5ty| 2 . Let us first proceed algebraically, considering possible perturbations to the Dirac equation which will 
open a gap, and only then ask how to obtain these perturbations from Stij. 

A given perturbation Stij will perturb the Dirac equation to 

Eip = (vf(a x k x + OLyky) + M^jip (5) 

where M is some matrix, the projection of Stij onto the space of the four states at the Dirac point. It may be shown 
that the perturbation M will be most efficient at opening a gap when it anti-commutes with a x ,a v . By efficient, we 
mean that we wish to maximize the gap for given Tr(M 2 ), as a first step to maximizing the gap for given |<5iy| 2 . 
Since there is only a sixteen dimensional space of matrices M, we can easily characterize all matrices that have the 
needed anti-commutation property; it is a four dimensional vector space. 

We will write three of the perturbations as matrices Mi, for i=l,2,3. In terms of a-matrices, they will be Mi = 
& Z ,M2 = f3,Ms = (3a x a y a z . These perturbations anti-commute with each other; in fact, one can make a change 
in basis in the Dirac equation which leaves a x ,a y unchanged, but produces continuous 0(3) rotations in the space 
of Mi, M2, M3. This continuous symmetry is only valid at low energy; it will be broken to a discrete symmetry by 
lattice effects as discussed below. By taking M = ^2 i rriiMi, for some numbers mi, mi, 771.3, we open a gap equal to 

V^2i( m i )■ We will refer to these as nonchiral mass terms. 

The fourth perturbation is of a different sort. It is M = m c M c with M c = ia x a y (here, c stands for chiral and we 
will refer to this as a chiral mass term) . This perturbation breaks parity and time-reversal symmetry. M c commutes 
with Mi, Mi, M3. As mentioned above, we are interested in the most efficient way for the system to open a gap; since 
M c does not anti-commute with Mi, it is most efficient for the system to take either purely chiral mass or purely 
nonchiral mass, so that m, = or m c = 0. 

Next, we would like to ask what perturbations in the £y will produce the desired mass matrix M. We will find 
that to produce Mi requires dimerizing the system by making the magnitudes of the \tij\ non-uniform; to produce 
M c requires adding additional fluxes to the system. Clearly, there is a large degeneracy here, as there is a only 16 
dimensional space of matrices M while in a given unit cell there is a 48 dimensional space of perturbations to tij. 
While some of the degeneracy is due to the large number of possible gauge transforms on tij , this does not completely 
alleviate the problem. Again, the question of efficiency becomes important: for each M, there is a class of tij which 
produce the desired perturbation, but only one element in the class minimizes Y) { ^ \tij\ 2 . We will find one unique 
perturbation in tij (up to arbitrariness in gauge) which produces the desired mass matrix. 
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The perturbation that we will pick for Ml is shown in figure 3. One can see that all the horizontal bonds have 
been either decreased or increased in strength, such that along a horizontal line the bonds alternate in strength while 
horizontal bonds which are in a vertical column all have the same strength. 

This perturbation is a dimerization of the iy . The spins form singlets most strongly across the largest ty , so that 
dimerization of iy tends to produce a spin solid and changes the long-range RVB to a short range set of singlets. 

This mass term breaks rotational and translational symmetry of the lattice. We will pick M2 and M3 to be lattice 
rotations of M\ . The continuous symmetry of the Dirac equation will then be broken at the lattice level to a discrete 
symmetry of permutations of mi, 7712, m 3 under lattice rotation, while lattice translations change the sign of any 2 of 
the 3 mass terms mi, 7712, m.3. We will find later that while we obtain symmetry breaking and produce a mass, the 
discrete nature of the lattice group will leave us with only pseudo-Goldstone modes. 

There are two mass terms which are symmetric under rotations. They are 

Mi + M 2 + M 3 

M i2 = ^ (6) 

-Mi - M 2 - M 3 , . 

Ms = -M i2 = l — 3 - (7) 

where 12 denotes the fact that the ty are strongest on the 12-site loop surrounding the unit cell, while 6 denotes the 
fact that the bonds are strongest on the hexagons and triangles. We show the perturbation to ty to produce M12 in 
figure 4. 

To estimate dimerization later, it will be useful to know connect the change in t to the eigenvalues of the mass 
matrix that arises. One finds that if the bonds on the 12-site loop are increased by 5t, while those on the hexagons 
and triangles are decreased by the same amount, then one produces a term mi2Mi 2 in the Dirac equation with 
mi2 = 1.57735....£t. 

For calculations later, it will be convenient to transform to a basis of gamma matrices. Defining j t — (3, 7« = Pan, 
we find that Mi = 73, M2 = 1, Ms = 275, while M c = i^tlxly 

It is interesting to compare to above characterization of possible perturbations in terms of gamma matrices to the 
situation in the 7r-flux phase on the square lattice, where there is again a Dirac spectrum, and again various mass 
terms can be introduced [2q |. We find that we are able to introduce one nonchiral mass term by dimerizing the 
horizontal bonds of the square lattice, so that the horizontal bonds alternate in strength as one moves horizontally 
along the lattice; another mass term can be introduced by dimerizing the vertical bonds. In the limit of extreme 
dimerization, these states correspond to short-range RVB states in which the dimers are stacked on top of each other, 
and all lie either horizontally or vertically 2(|. By taking sums of these two mass terms, we can produce a short 
range state in which dimers resonate around a square. The final nonchiral mass term can be obtained by placing an 
on-site potential on one sublattice of the square lattice; this corresponds to introducing Neel order into the system 
pj| . Due to the highly frustrated nature of the kagome lattice, in this paper we will not have any such terms involving 
introducing on-site potentials. 

For the square lattice, the chiral mass term can also be introduced. It requires adding additional couplings to the 
system, which connect diagonally across a given plaquette, and then inserting 7r/2 flux through the triangle that is 
formed when a particle traverses two sides of a plaquette and then crosses the plaquette on the diagonal |2q]. 



III. CONNECTION TO OTHER VALENCE BOND STATES 



From the parent state, it is possible to continuously connect to other possible valence bond states, using the mass 
terms we have found above. Let us first consider the chiral spin liquid and then the spin-Peierls solid of Marston and 
Ze n g||. 

Let us consider a state such that the flux through each triangle is equal to 8, and the flux through each hexagon is 
equal to 7r — 28. Then, as 8 varies from to \ we continuously deform from the parent state to the chiral spin liquid. 
By looking at how ty changes along this deformation, and then projecting this change onto the Dirac point, we find 
that for small 8, the perturbation exactly produces the chiral mass term M c ! Let us note that at 8 = j there is a 
highly interesting band structure, discussed in the Appendix, with multiple flat bands. 

The chiral spin liquid state improves on our parent state at the mean-field level. Since we will argue below that 
fluctuations stabilize our state against the chiral mass term, let us here analyze why the chiral state works at the 
mean-field level, and provide a qualitative explanation of why fluctuations destroy the chiral spin liquid. 

The idea behind the chiral state results from the "Rokshar rules" |29), which argue that one should put flux ^ 
through every triangle, no flux through the hexagons, and have a total flux of it through the loop of length twelve 
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that surrounds a hexagon and its six attached triangles. These rules are derived from considering individual hexagons 
and triangles in isolation, and minimizing the mean-field energy. 

While our parent state appears to violate every one of the rules, except the rule for the length-12 loop, the chiral 
spin liquid is in perfect agreement with these rules! Let us focus on one isolated triangle, with \t\ = 1. If there is 
no flux through the triangle, there are two negative energy states with energy —1. We can put two particles in one 
state, and one in another, for a total energy of —3. By adding ? flux, we have one state at energy — v3 and another 
state at energy zero. By putting two particles in the lowest energy state, we improve the energy of the system, and 
introduce a chirality. 

Now, consider the triangle coupled to the rest of the system. If the rest of the system strongly scatters the particles 
in the triangle, it may no longer by appropriate to think of two particles in one state and one in another. One 
must instead think of each of the two negative energy states of the triangle as each having average occupation of 
one-and-a-half particles. In that case, it is most advantageous to put no flux through the triangle. 

So, if the system to which the triangle is coupled is chiral, so that all triangles in the system have the same flux 5, 
then the chiral spin liquid may work. But if the triangle is coupled to states which are not chiral, then the chiral spin 
liquid is destroyed. One sees this even at the mean field level; a state in which triangles have alternating flux ±f is 
significantly worse in energy than our parent state. Similarly, if one introduces sufficiently strong dimerization mi2, 
one finds that the system is stable against weak m c . Within the RG below, we will consider the fluctuations in the 
masses rrii and show that they help stabilize the parent state against the chiral perturbation. 

The spin solid of Marston and Zeng can also be obtained from the parent state. In this spin solid state, the idea is 
to look for dimer coverings which maximize the number of "perfect hexagons" , hexagons on which three of the bonds 
are covered by dimers. Attached to these hexagons are "defect triangles" , triangles on which no bonds are covered by 
triangles. Clearly, we wish to increase \tij\ on the perfect hexagons, while decreasing it on the defect triangles. This 
will project onto the mass term Mq on the 12-site cell that includes the given perfect hexagon. 

The perfect hexagons are then supposed to form a lattice. We can obtain this lattice by taking the triangular 
lattice of 12-site cells, and placing perfect hexagons inside the 12-site cells on two out of the three sublattices of the 
triangular lattice. In this case the cells containing perfect hexagons form a honeycomb lattice. This gives rise to a 
staggered mass state. We produce a mass term Mq on two-thirds of the system, so that the Dirac particles feel a net 
Mq at zero momentum, as well as a fluctuating Mq at finite momentum. 

IV. NAIVE MEAN-FIELD AND PROJECTED MEAN-FIELD 

To compare the energies of possible RVB states, including the various mass perturbations of our parent state, we 
turn to the RVB mean field. When we look for a mean- field solution of in equation (||), it is known |2l| that the 
best mean-field is a dimer covering. Still, let us start by looking at results of the naive mean-field calculation, and 
then later provide a projected mean field calculation. 

Let us introduce the Green's function between sites, Gij, defined to be the sum over all occupied fermionic states, 

-0, Of 

(8) 

The fermionic energy is then equal to 

2 (-'.j 1 - (9) 

<i,j> 

where the factor of 2 arises from the presence of up and down species of fermion. 
For our parent state, explicit calculation shows that, for nearest neighbors 

\Gij\ = 0.221383... (10) 

From the mean-field condition for equation (@), we find for the parent state that t = 0.221383 J, so by equation (0) 

v f = (0.0904.. )J (11) 

We find that within the projected mean-field that the parent stable is unstable to all of the massive fluctuations, 
including the chiral mass fluctuation which will drive the system to a chiral spin liquid. 

For infinitesimal perturbations, the different nonchiral masses all provide an equivalent improvement in mean- 
field energy. To some extent this is due to the approximate low energy symmetry of the Dirac equation to rotating 
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continuously among the different mass terms. However, it is interesting that lattice effects do not break this symmetry. 
The reason is the discrete lattice symmetry. The change in energy for taking M — m\M\ + m^Mi + 1713M3 is, for 
small m, a quadratic form in m,-. Let this form be 



E 



Cijrrnmj (12) 



The coefficients 011,022,033 in this form must all be the same due to lattice rotation symmetry. Lattice translation 
symmetry permits one to change the sign of any two of the m%, and prevents a nonvanishing cy- for i =/= j. Therefore, 
for small perturbations the energy gain for introducing a mass must be dependent only on the magnitude of the mass, 
not the particular mass term used. For larger perturbations, the energy gains may depend on the particular mass 
term used, and of all the nonchiral mass terms, the system gains the most energy by a mass M\%. 

Now let us turn to the projected mean-field. Instead of doing a full Gutzwiller projection, we will use an approxima- 
tion introduced by Hsu [ p7| . Within a variational Gutzwiller projection, one minimizes the energy of the Hamiltonian 
(|l|) . Hsu's idea at the lowest level of approximation is to note that the Hamiltonian is a sum of terms J Si ' Sj over 
different neighbors i,j, and, when evaluating the expectation value of each of these terms, to perform the projection 
only on the given sites At this level, the variational principle corresponds to minimizing 

<hj> <i,j> 

over all possible tij , where Gy is determined by ty . 

For our parent state, we find that ~6 1+16 ]g. |4 = —0.2832... By going to the chiral spin liquid, the system improves 
the ground state energy by roughly 2.9 percent within the projected mean-field approximation. By going to a state 
with staggered ±? flux through each triangle, the systems worsens the ground state energy by roughly 1.4 percent. 
Within this approximation the system is stable against the nonchiral mass perturbations as all the nonchiral mass 
perturbations worsen the ground state energy at this level of approximation. Again due to discrete lattice symmetry, 
the energy cost is independent of the particular mass term for small mass, while for larger perturbations, the energy 
costs differ, and the M12 perturbation costs the least energy. 



V. FIELD THEORY OF FLUCTUATIONS 



Having discussed the naive and projected mean-fields for the problem, we must include fluctuations about the 
mean-field. To do this, we will use an SU(N) generalization of the original problem j30), such that the projection 
procedure of Hsu becomes exact. After discussing how to do this in the abstract, we will present the field theory for 
our specific problem: it will have a number of interacting fields, including fermions, gauge fields, the nonchiral mass 
terms discussed above, which we will refer to as "pion" fields, and the chiral mass term. 

The approximation of Hsu ammounts to minimizing equation (|l3|). By introducing an auxiliary field A^ we can 
"decouple" this sum of functions of Gij and instead extremize the function 

j2 G i^i+m»\ 2 ) (14) 

where /(|A| 2 ) is a Legendre transform of pjjjgWi ■ 

Then, we can interchange the order of extremizations, and extremize this quantity over ty before extremizing over 
Xij. We find that this is extremized at iy = Xij. Then we proceed to extremizing over the one remaining set of 
variables Ay . But since tij = Xij , we are are equivalently trying to extremize the function 

E G «*ii+/(fe'| 2 ) ( 15 ) 

over all fy. Note, now, that Gijtji is exactly the kinetic energy of the fermions. So, finally, we are trying to extremize 

h = ]T (rimjMj) + h.c) + ]T /(|M 2 ) (is) 

<i,j> <i,j> 



G 



Returning to the language of functional integrals, we can introduce an SU(N) field theory for which the Hsu 
projection procedure becomes exact. We take a large N limit in the number of fermion fields in equation ([If]), so that 
a = 1...N. Then, we integrate over all possible iy in that equation, undoing the decoupling procedure above, and 
rewrite the result in terms of spin operators. We find 

h= y Si ' s > * a?) 

We should note a few facts about this procedure. When we demonstrate the equivalence of the large N mean- 
field with the Hsu mean-field, it is the large N limit that permits us to ignore fluctuations in fcy,Ay, so that the 
decouplings amounts exactly to taking a Legendre transforms; at finite N, the decoupling of an interaction is not 
exactly a Legendre transform. Further, we used the word "extremize" above with care: in some cases we maximize 

IGI 2 

while in other cases we minimize, as in some places the function 1+16 ' G |4 has positive curvature while in other cases 
it has negative curvature. This does not provide any formal problems when performing the decoupling at the level of 
functional integrals, so long as we correctly choose the integration contour of Ay. 

The fractional operator, X)<« j> in equation ( |l7|) may be interpreted as a formal power series, so that it 

includes operators of the form (Si ■ Sj) k for all k. At N = 2, this operator is equivalent to Si ■ Sj, up to a constant 
factor. 

Finally, the above procedure is similar to the technique of introducing biquadratic interactions, (Si ■ Sj) 2 , into the 
Hamiltonian to stabilize RVB states against dimerization. We simply prefer the above Hamiltonian as it reproduces 
exactly our desired mean-field theory. 

Having defined a large N theory with no fluctuations in the decoupling fields, we next add in fluctuations. Formally, 
this can be handled by a 1/N expansion. We will directly write the field theory at TV = 1 without including explicit 
factors of N. 

The theory includes several modes. There is the Dirac fermion field, ip(x). This is coupled to a fluctuating U(l) 
gauge field, A^(x), fi — t,x,y. By considering other fluctuations in ty we will also obtain a fluctuating chiral mass 
field that we will refer to as u(x), and a triplet of fluctuating nonchiral mass terms that we will group into one 
"pion" field -K a (x), a = 1, 2, 3. If the pion field acquires an expectation value, then the fermions will acquire a mass 
m a = (ir a ). 

For the field theory, we will suppress the velocity v/. At the level of the bare action the pi, sigma, and gauge fields 
can have different velocities from the fermi fields. However, the greatest contribution to the effective action of the 
bosonic fields arises from integrating over the relativistic fermions, so that at low energies the velocity of the bosonic 
fields must be roughly equal to that of the fermionic fields. 

The Lagrangian L we will take is 

L = J d 3 xL f + L A + L M (18) 

La = TTa V" ( 2 °) 
JL^a(x)(dl + ml)a(x) + J-^)^ + (m^ 6 )^) (21) 



where 



J->M 

We have written the mass for the pion field as a matrix. However, following the arguments for equation ([L2|), the 
masses of the different pion modes must be the same. It is only after condensation of the pion field, breaking the 
lattice symmetry, that the masses can differ. 

We have inserted factors of A, representing a lattice cutoff scale, into the action to make the coupling constants 
dimensionless. We have chosen to scale the bosonic fields so that all coupling constant dependence appears in the 
action La and Lm, not Lf. 

In addition to the terms we have written, there must be a quartic interaction term for the n and a fields. This term 
is necessary to stabilize the action if the system spontaneously breaks a symmetry and has either < or m ff < 0. 
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Since we will be initially starting the renormalization group of the next section with both such masses positive, we 
can temporarily ignore the quartic term at high energy under the assumption that this term is small. If the system 
acquires an expectation value for the it fields, the quartic term will break the continuous symmetry down to the lattice 
symmetry, and give a small mass for the approximate Goldstone modes. There will also be cubic terms that give 
a mass to these modes. 

Another interesting term we have left out is ^j^^F^ip, which can be added to change the g-factor of the Dirac 
fermions. In the absence of external fields, there are two degenerate states of the Dirac equation at each energy. For 
physical electrons, this reflects a spin degeneracy. For the spinons we consider, which already have a definite spin, this 
degeneracy instead reflects a chirality degeneracy, and we will refer to it as such. If the system has an odd number of 
sites, and hence an unpaired spinon, not only does the system have a net spinon, but it also has a net chirality, which 
can be taken to be positive or negative. Generically the <?-factor will be non-zero. 



VI. RENORMALIZATION GROUP 



We will consider the one-loop RG from the field theory. We will see that it is indeed possible for fluctuations to 
lead to a condensation of the pion field. 

For the gauge, pion, and sigma fields we will use a simple mode elimination RG, with a cutoff A. For the fermion 
fields, we will introduce a set of massive regulator fields with masses of order A and reduce the regulator mass. The 
choice of the particular mass terms for the regulator fields represents a lattice breaking of the Goldstone symmetry. 
It is possible to preserve the needed lattice symmetry of equation ( |l2| ) by introducing seven regulator fields. Four are 
ghost fields with masses proportional to Ml + M 2 + M 3 , Mi - M 2 - M 3 , -Mi + M 2 - M 3 , -Mi - M 2 + M 3 and the 
other three are not ghosts and have masses proportional to M\,M%^ M3. 

Initially the theory will have a cutoff Ao, defining the lattice scale. As we renormalize, we lower the cutoff A, and 
rescale all distances and fields to keep A fixed. 

We must take into account self-energy corrections to the fermions from interactions with the bosonic fields, vertex 
corrections, and self-energy corrections to the bosonic fields from vacuum polarization bubbles. If we were to take 
into account only the self-energy corrections to the bosonic fields, and not the vertex and fermionic self-energy terms, 
we would find that we are considering just the mean-field theory in the bosonic fields. 

We find the following RG equations: 



d^9A = _ 2 9a_ 
cfln(A /A) 3 



(22) 



dlng n 3g\ -\-g\-g\ . 

= H 1- vacuum polarization (23) 



efln(A /A) 2tt 
dlng a 3g\-3g%-gl 



<fln(A /A) 2tt 2 
dml h 



4- vacuum polarization (24) 



dln(A /A) 
dln(A /A) 



2 + vacuum polarization (25) 



= 2 + vacuum polarization (26) 



We have avoided explicitly writing the vacuum polarization contributions to the sigma and pion fields. The vacuum 
polarization contribution to the mass is regularization dependent, while the vacuum polarization contribution to the 
coupling constant is ultraviolet convergent and is dominated by the infrared contribution. 

Fluctuations in the gauge field increase the coupling constants for the pion and sigma fields. This reflects the 
binding force due to the gauge field between charged spinons, and the resulting tendency to break chiral symmetry. 
Further, we see that the coupling constant for the pion field increases more rapidly than that for the sigma field, 
reflecting the destabilization of the chiral state by fluctuations in the pion field. 

Thus, we see from the renormalization group that there is a range of bare parameters such that the theory will 
condense the pion field, producing a nonchiral mass term for the fermions, even though at the mean-field level the 
theory would rather produce a nonchiral mass term for the fermions. 
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In the next section we will consider the low energy action after condensation. We will first discuss the mass term 
for the fermions that appears. 

Unfortunately, it is beyond our ability to calculate the bare parameters in the field theory with any precision, and 
so the mass of the fermion field is not something we can compute. Let us instead take the mass of the fermion as 
an input from numerics, and use that to check for consistency of our theory. Extrapolating finite size results from 
systems of up to 36 sites, one finds that the system has a gap to triplet excitations which is of order J/20 or less [p4[ . 
While the gap is decreasing even at the largest sizes, it appears to be bounded below by roughly J/40. Assuming 
that the triplet excitations are made of pairs of spinons, we can estimate the spinon gap as being half the triplet gap. 
Further evidence for the spinon gap being roughly half the triplet gap comes from odd-even studies of the energy 
dependence on N p4fl . The fermion mass is half the spinon gap, or one quarter the triplet gap. 

Using this estimated spinon gap, and the calculated velocity of our Dirac particles from RVB theory, we can 
obtain the correlation length of the Dirac particles. Taking the estimate of J/20 for the triplet gap, we find that the 
correlation length is roughly 8 of our twelve-site unit cells, large enough to include the N — 36 numerical studies. 

We can also estimate the strength of dimerization at the mean field level, by asking how large a change in Uj 
is needed to produce the desired mass term. Assuming that the dimerization is provided by a perturbation of the 
form M12, one finds that the £y 12-site loops are increased by approximately 3.5 percent, while the other bonds are 
decreased by 3.5 percent. This is a relatively small amount of dimerization, and we expect that only after a significant 
increase in system size will numerical studies be able to detect this directly from a dimer-dimer correlation function. 



VII. LOW ENERGY MODES 



The remaining low energy modes after the pion field condenses are the Goldstone excitations of the pion field and 
the gauge excitations, which we will argue provide the low energy singlet modes seen in numerics. While numerical 
calculations have only probed systems up to N = 36 sites, which is relatively small considering that we take a unit 
cell of 12 sites, experiment also reveals a quadratic low temperature specific heat. This specific heat suggests that a 
bosonic mode with a linear density of states survives to much larger scales, while the insensitivity of the specific heat 
to weak magnetic fields suggests that the mode is still a singlet. In this section we will first discuss the nature of the 
low energy modes and then the ultimate fate of our pion and gauge excitations at large distances, including a gap for 
the pion from lattice effects as well as a confining phase for the gauge fields. 

Once the pion field condenses, the system is left with two pseudo-Goldstone pion modes as well as gauge modes. 
The gauge action is 



1 

4Agj 



L a = TTjtV" (27) 



with 



91 oc m (28) 



where m = \{tt)\ is the fermion mass. The pion action will be a sigma model. If we change the normalization on the 
pion field so that |tt| = 1, we get the model 



L w = ^(M 2 + ^ttW ,2!,, 



where the coupling constant g 2 is proportional to to -1 and the mass m^r represents the breaking of continuous 
symmetry by lattice effects. 

We have chosen the mass term for the pion to cause the pion to prefer to condense in a way that gives rise to a mass 
M\2- In the projected mean- field calculation above, we considered states invariant under the rotational symmetry, so 
that \mi\ = | Trial = l TO 3|- This provided two inequivalent perturbations. In one, we increased \Uj\ on hexagons and 
triangles; in the other we increased \tij\ on a loop of length twelve. While at the mean-field level symmetry breaking 
does not occur, from the projected mean-field calculation we can still argue that the preferred symmetry breaking 
pattern would be given by a mass matrix M±2- Other patterns are of course possible, and comparison with numerics 
provides some evidence that a staggered mass is also a possibility, at least for small systems; in the conclusion we 
discuss possibilities of numerically testing the preferred mass pattern. 

In the continuum field theory, the pion mass (m 2 r ) ab is ultraviolet divergent. However, lattice symmetry forces the 
masses of the pion modes to be the same before condensation. In order to use the continuum theory to estimate 
after condensation, we need to turn to the cubic interaction terms in n. These are 







d X537ri7r 2 7r 3 (30) 

with a cubic coupling constant <?3 that is generically of order unity. Inserting an expectation value of it of order m, 
we obtain a quadratic term in 7r. Including this quadratic term in (mj) ab , this will cause the masses of the different 
pion modes to differ by order m so that m^r will be of order m. 

However, we can obtain a better estimate of the mass difference numerically from the the projected mean-field 
calculation of energy; while we did not obtain pion condensation at the mean-field level, for a given expectation value 
of the pion field, we can use the difference in mean-field energies for various continuous rotations of the pion field to 
obtain an estimate of the pion gap. Using the numerical estimate for the triplet gap, and hence the estimate for the 
fermionic mass m, we have calculated the projected mean-field energies for taking M = mM\i and M = vtlMq. The 
difference in energies is 0.000396 J per 12-site cell, so that m^r w V0. 000396m J. This is small enough that we can 
ignore this mass for most purposes. Evidently, the cubic coupling constant g 3 is very small. 

While the pion is already gapped by lattice effects, instantons will gap the gauge field, leading to confinement of 
the spinons. The gauge field describes compact QED in 2+1 dimensions, which is confining for all gA |3l[]- The gauge 
coupling is proportional to m, so that the action for an instanton will be of order 

Soc^ (31) 
m 

The instanton density is proportional to e~ s . In the weak coupling limit, the instantons lead to a gap for the gauge 
field of order Aqe~ s . As this is exponentially small, we can ignore the gap in the gauge field. 



VIII. FINITE SIZE SYSTEMS AND CHERN NUMBERS 



In this section we will consider some effects of finite size systems to begin comparison with numerics. First we will 
consider some complications in defining the parent state on systems with an odd number of sites, which force the 
system to have some net flux. Then we will show how this leads to a degeneracy in the spectrum and nonvanishing 
Chern numbers for the states under twist in the spin boundary conditions. Finally, we will discuss some effects of 
finite size for even size systems. 

One of the most interesting results found in numerical studies of odd size systems is a non- vanishing Chern number 
pi| for the ground state of ±1. This is a quantity that provides an analogue for a spin system of the quantum Hall 
effect p2j . Since the Hamiltonian of equation (Jl]) does not explicitly break time- reversal symmetry, a non- vanishing 
Chern number requires a spontaneous breaking of time-reversal symmetry. However, the spontaneous breaking of 
time-reversal symmetry is not enough, as other spin systems that break this symmetry have vanishing Chern number 
J32| ; the kagome antiferromagnet may be the first Hamiltonian with time-reversal and parity symmetry where a non- 
vanishing Chern number has been observed. In order to understand the appearance of the Chern number, we must 
first understand how to form the parent state on an odd size system. 

The systems studied numerically have periodic boundary conditions, so that they live on a torus. On system defined 
on a torus the net flux penetrating the surface must be an integer multiple of 2tt. One can also add solenoid fluxes, 
61,02, defining the phase that the spinon acquires when traversing a topologically nontrivial closed loop around the 
torus. For simplicity, we will use coordinates on the torus which range from to 2tt in both directions, although in 
actuality for the kagome lattice the systems considered are not square. 

The parent state has it flux through each hexagon. On a system with N sites, there are N/3 hexagons, and so on 
a system with an odd N, one would like to have a net flux through the system that is an odd multiple of n. This is 
not possible, and so the system must have some additional flux so that the total is a multiple of 2-7r. For example, on 
a system with N = 27, there are 9 hexagons, so the system can put ir ± J flux through each so that the total flux is 
either 87r or 107r. The system then must become chiral and break time reversal symmetry since it cannot construct 
the parent state. 

A qualitative way of describing this effect is to say that for a system with an odd number of sites, there is an 
unpaired spinon, which has a chirality. The spinon then couples to the gauge field and produces a flux. 

Given this net flux through the system, let us consider the effective Dirac equation for the spinons. The results we 
get for the Chern number do not rely on the Dirac description, and can be derived directly from the lattice model; 
we feel that the Dirac method is more elegant and gives more physical insight. 

The Dirac particles feel the extra flux that has been added, and so the spinons move in a magnetic field, such that 
the net flux the spinons feel is exactly ±7r. Again, there seems to be a contradiction, since it is not possible for the 
system to have a net flux of ±7r flux through the torus. The answer to the contradiction is that the Dirac particles 
have an extra chirality index. So, in addition to including solenoid fluxes for the Dirac equation, the Dirac particles 
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can change chirality when completing a loop around the system. Let us then generalize the solenoid flux to a pair of 
4-by-4 matrices t/i, U2, describing the change in the wavefunction when the particle completes a loop. 

Then, when the particle traverses a loop around the torus from (0,0) to (0, 2tt) to (2n, 2ir) to (2tt, 0) to (0,0), the 
wavefunction gets multiplied by 

- UxU 2 Ulul (32) 

where the minus sign is from the magnetic flux through the torus. Since equation (^) must be equal to 1, we find 
that U1JJ2 necessarily anti-commute. 

One may regard the matrices U\ , U2 as arising from a non- Abelian gauge field connecting opposite chiralities of the 
spinons. The commutator of the matrices represents an additional flux of it from the non- Abelian field, giving a total 
flux of 2tt on the torus. The extra n flux from the non- Abelian field is the flux that arises from having an odd number 
of hexagons on the lattice, so that when the particle completes the given loop around the lattice it has enclosed an 
odd number of ir fluxes. One sees that the non- Abelian flux is localized at a point, although one must be careful that 
this localization at a point does not imply a breaking of translational symmetry. 

The addition of matrices U±,U2 is natural from the lattice point of view. The unit cell which includes both chiralities 
of Dirac particles is 12 sites, while the smallest unit cell possible for the parent state is 6 sites. Since there is no 
way to cover an odd size lattice with 6 or 12 site unit cells, something must scatter between chiralities, as when the 
particle completes a loop it has changed between chiralities. 

To give a very simple example of this, consider a one-dimensional ring with an odd number of sites. The natural 
unit cell for the a one-dimensional chain is two-sites, to include both Dirac points. If the particle moves around an 
odd-length ring, two sites at a time, it must return to the starting point displaced by half a unit cell. 

To give a slightly more complicated example, consider the 7r-flux phase of the square lattice (33) for a system of 9 
sites, shown in figure 5. Solid lines represent bonds within the cell of 9 sites, while dotted lines represent bonds to 
provide toroidal boundary conditions. There are 9 squares in the system, and so there will be 7r± ^ flux through each 
square. 4 of the squares lie within the cell and are labeled 1-4, another 4 lie to the sides and are labeled 5-8, while 
the 9th square lies in the corner. The natural unit cell for the Dirac particles is 4 sites, so when a particle completes 
a loop around the torus it is displaced by half a unit cell. The wavefunction is multiplied by a matrix U\ for a loop 
in the i-direction and a matrix U2 for a loop in the y-direction. Precisely due to the odd number of squares, one 
finds that the matrices U\, U2 anti-commute. It is natural to think of the non- Abelian flux as arising from the 7r-flux 
through the 9th square, on the corners. 

Returning to the kagome lattice, let us now look at the wavefunctions of the Dirac equation with this magnetic 
flux. It is convenient to find the wavefunctions by enlarging the torus by a factor of two in each direction, as shown 
in figure 5. The + and — symbols in the figure denote the chirality of the particle in each quadrant. When the 
particle completes a circuit on the original torus, it changes chirality, and hence moves into a different quadrant of 
the enlarged torus, while picking up a phase. 

The net flux on the enlarged torus is equal to 4 times the flux on the original torus, or An. One might imagine that 
there will be extra sources of 7r-flux on the enlarged torus at the points where the quadrants meet. However, since 
the 7r non- Abelian flux is purely a result of an odd number of hexagons on the original torus, we can drop the extra 
sources of flux on the enlarged torus, and we are left with an explicitly translationally invariant problem of a Dirac 
particle moving in a constant magnetic field. 

On the enlarged torus, the Dirac equation becomes a two-component equation 

(iEa z + a x (A x + id x ) + a y (A y + id y )^ijj = (33) 

Taking the square we find 

E 2 4> = ({id x - A x f + (id y - A y f + a z B)4> (34) 

This is the well known equation for a Dirac particle in a magnetic field, and has Landau levels. 

Since the net flux through the enlarged torus is equal to 47r, there are two-states in each Landau level for each tr z , 
hence four states for each Landau level in total. As we are dealing with a two component equation, only one sign of 
E is allowed in equation ([54j) for given a z . 

Therefore, the energy levels on the enlarged torus are doubly degenerate. However, the enlarged torus has an 
unphysical degree of freedom: opposite quadrants describe the state of the particle on the original torus. So, the 
energy levels on the original torus are only singly degenerate and the spectrum is discrete with one level at zero 
energy. This is the relativistic generalization of Landau levels. For an odd size system, all Landau levels below the 
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zero energy are occupied, and hence filled, for both spin up and spin down particles, while the zero energy level is 
occupied only by one unpaired spinon. 

Numerically, the ground state of the many-body system has been seen to have an extra degeneracy factor of two, 
beyond the trivial spin degeneracy. This is a consequence of the spontaneous generation of the magnetic field, so that 
the system can pick either sign for the field. 

An interesting way of viewing the Landau levels is that we have ir Abelian flux, implying that there are 1 /2 states 
per Landau level. Multiplying the 1/2 by a factor of two for chirality degeneracy, we get one state per Landau level. 

We can now introduce the Chern number of the system, which characterizes a transverse response of spin currents. 
Let us adjust the boundary conditions of the system so that 

S ± (x,y)=e ±i ^S ± (x + 2ir,y)) (35) 
S ± (x,y) = e ± ^S ± (x,y + 2ir) (36) 

where 4>i,4>2 are angles. 

If the ground state is a wavcfunction ^ , then the Chern number is defined as the integral 




This number is quantized, and non-vanishing only for complex states. Since the Hamiltonian does not break time- 
reversal symmetry, complex conjugate states are degenerate with opposite Chern numbers. 

In the presence of these boundary conditions, the spinon boundary conditions, with additional self-generated fluxes 
9i,02, become 

1>l(x,y)=j$il,i(x + 2Tr,y) (38) 



ipl(x,y)=e i ^4(x + 2ir,y) (39) 

where 

p* d = (20i±<h) (40) 

and similarly for the other direction. The 2ir periodicity in <p> is not obvious from equation ( |38| , |3"9| ) , but the ability of 
the system to adjust 9 produces the desired periodicity. 

To give a simple example of how a system can adjust 8, consider a system of 4 sites on a ring. In the absence of 
a twist in boundary conditions 0, the system places 9 = ir flux through the ring. As <f> increases, 9 remains equal 
to 7r, and the mean-field energy of the system gradually increases until <p = ir. At this point, 9 jumps to 0, and the 
mean-field energy begins to decrease for increasing <p>. So, as <p varies from to 2ir and 9 jumps as described above, we 
find that | = ( 2fl +^) varies from ir to ^ to ^ to ir. In order for 9 to jump like this, the spinon states with | = ir ± ^ 
must be degenerate. 

Now, we can look at the Chern number of the system, assuming non-interacting spinons. It then amounts to a 
Chern number calculation of the fermionic wavefunctions. Assuming non- interacting spinons, we can get the change 
in $ in equation ( |37|) from the change in the spinon wavefunctions. 

While in general the wavefunction gets multiplied by a matrix on moving around the original torus, only the U{1) 
part of this matrix adjusts in response to changes in <p>. The U(l) part of the matrix is just the angle 9. Carrying out 
the calculation on the enlarged torus, we find that the boundary conditions become 



4>Ux,y) = e lp m(x + 4ir,y) (41) 



and similarly for down spinons. 
So, we wish to compute 



2W J "d P ^ d ydp 2 l ' d y 

summed over all spinon wavefunctions ip. 



1 f f dip dip dp"' dp\ 



(42) 
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In equation ( pl| ) the periodicity in <fi seems obvious even without 6, as on the enlarged torus the periodicity of 
the spinon wavefunctions in response to a twist in boundary conditions is halved in both directions. However, we 
have introduced the degeneracy of two on the enlarged torus representing the fact that on the original torus the 
wavefunctions are periodic in (pi, P2) with periods (0, 2ir) and (jr, n) but not with period (0, tt), and as a result only 
half the possible wavefunctions on the enlarged torus are physical. 

In order to keep the wavefunction in the physical sector, 9 must jump discontinuously by n as <f> changes, and 
as a result for a given spinon state we only integrate equation ( fl2] ) over half the torus of possible phases (pi,p2)- 
The fact of integrating over half the torus, or equivalently the fact that the Landau levels contain one physical and 
one unphysical wavefunction, does not prevent a defining of the Chern number for the spinon wavefunctions. When 
6 jumps it connects two degenerate states, and so the contribution of equation ([l2|) to equation (^) must still be 
quantized as an integer which we can still refer to as a Chern number for the spinon. 

Within the lattice formulation there are no further conceptual problems and we must simply compute the integrals, 
but within the continuum Dirac equation we must account for the negative energy sea. The correct understanding of 
this was found by Haldane Qj. 

One must add a massive regulator field, and compare the difference in Chern number between the massless and 
massive fields. The massive field has the same Landau level spectrum, but no zero mode. So, the difference in Chern 
numbers is due to the zero mode, which sits in the lowest Landau level. There are two states in the lowest Landau 
level, one physical and one unphysical. It is the physical state that carries the Chern number of ±1, giving the ground 
state of the spin system a net Chern number of ±1, as seen numerically 0]. 

We expect that the low-lying states will continue to have an odd Chern number, in agreement with numerical 
results. If a particle-hole pair is excited within the Dirac band near the Dirac point, the Chern number will not 
change. If the particle is excited from the band edge, the Chern number can change by ±2. Only if a particle is 
excited from the flat band to the Dirac band can the Chern number change by ±1, giving rise to an even Chern 
number. However, these states will be much higher in energy. One can also consider excited states with net Abelian 
flux be equal to 37r, 57r, ... 

It is very interesting to think about these possible excited states which may have more than 7r flux for the Dirac 
particles. The two-component particles carry a a z index, which will couple to the magnetic field. If a large field is 
induced, a number of spinons of the same a z will be produced in the zero mode, so that the total number of spinous 
in the zero mode is odd. For one spinon we had one filled Landau level, with one particle. With several spinons one 
might be able to construct fractional Hall states of spinons. 

We have argued that in the thermodynamic limit the system will acquire a mass. On an odd size lattice, the mass 
term must change sign somewhere, as the lattice cannot be tiled with 12-site unit cells. At the domain wall where 
the mass changes sign, one expects to trap a midgap state, so there still should be a zero mode, even with mass. This 
may permit the nonvanishing Chern number to survive. 

Ret urning to even system size, let us consider the size dependence of the triplet gap. The energy of the spinon is 
E = yj (vfk) 2 + to 2 . In the absence of a solenoid flux, the smallest k would be equal to zero, but by creating a solenoid 
flux the energy can be improved and the smallest k will be of order the inverse linear dimension of the system, or 
N^ 1 / 2 . As a result, the triplet gap is decreasing with system size, in agreement with numerics. By twisting the spin 
boundary conditions one may be able to reduce the gap to S z = ±1 excitations. It would be interesting to look for 
this effect. 

Further, in the presence of these solenoid fluxes, other fermionic states with k ^ will become approximately 
degenerate with the k — fermionic state. This means that a spatially varying mass term which scatters between k 
states can open a gap just as well as the spatially constant mass term can. This will be important when we consider 
the low energy Goldstone excitations on finite size systems, below. 



IX. GOLDSTONE MODES, TOWER STATES, AND NUMERICS 

After breaking a symmetry, and giving mass to the fermions, the system is left with low energy pion and gauge 
modes. Above, we argued that the gap for these modes is too small to be seen in numerical calculations. In this 
section, we will treat these modes as gapless and discuss the energy spectrum that results for finite size systems to 
compare with numerical calculations. 

It is known ]35[ | that breaking a continuous symmetry gives rise to two kinds of low energy modes. First, there 
are the Goldstone modes with non-zero wavevector k. In the case of our pion and gauge modes, the energy is then 
proportional to k. For a 2 + 1 dimensional system with N sites, the lowest k is proportional to TV -1 / 2 and so the 
lowest Goldstone excitation has energy proportional to iV -1 / 2 . 
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Second, there is the "tower" of k — modes. These correspond to global rotations of the entire system, and have 
an energy proportional to N^ 1 . 

Numerical diagonalization [ [36| of the triangular lattice Heisenberg anti-ferromagnet, which has Neel order, shows 
very clearly the distinction between the tower of states, and the k ^ states (spin waves). However, no such distinction 
is found in the kagome lattice 0, no separation between low energy modes of energy TV -1 and TV -1 / 2 . Within our 
model, this is to be expected for N = 36. Since the effective action of the pion and gauge fields arises from integrating 
out the fermions, this action must be approximately relativistic, with the same velocity as the fermions. So, even 
without explicit calculation, we can obtain the energy of the lowest k ^ mode directly from the velocity appearing 
in Dirac equation. For the largest numerical diagonalizations, systems with N = 36 total sites or 3 of our 12-site unit 
cells, this energy turns out to be of order the triplet gap, and so this Goldstone mode is too high in energy to appear 
in the continuum of low energy singlets. 

The only states that will be observed in numerics are states in the tower. We can obtain the energy of these states 
from equation (p9|), assuming that tt(x) is constant. Then we get, assuming small m^, 

L=^f(d t n«) 2 (43) 
9 

where NA^ 2 is the area of the system. The states of this will be spherical harmonics, perturbed by mass term. With 
g 2 cx m , these states will have an energy proportional to 



Nm 



(44) 



and so for sufficiently large N will be below the triplet gap. A more precise knowledge of the prefactors will be needed 
to determine whether this is low enough to correspond to the low energy modes seen numerically. 

There will also be "tower" states for the gauge field, which correspond to different solenoid fluxes through the 
system. The energy for these can be obtained from equation ( p7| ) assuming that is constant over the sample. To 
get the energy for these we have to remember that the gauge group is compact, and realize that F^ v is derived from 
a set of {7(1 ) matrices with a lattice length A . Then, the energy of the gauge states is proportional to 

which is definitely below the triplet gap and certainly small enough to be the origin of some of the low energy modes 
in numerics. 

There is one puzzle involved in the tower states. It was observed numerically that, on 36 site samples, the energy 
of the lowest energy state of the system at given total momentum did not vary appreciably across the Brillouin zone 
p4j. This may seem to be in contradiction to the hypothesis that the low energy states come from the tower. The 
resolution of this may lie in realizing that there are only 4 inequivalent points in the Brillouin zone of the system. 
One of these is at k = 0, while including a constant nonvanishing mass term M\% also breaks translational symmetry 
on the kagome lattice, and can give one more point in the Brillouin zone. To obtain the last points in the Brillouin 
zone, one must remember that for small systems there can be other symmetry breaking patterns with staggered mass, 
breaking translational symmetry in different ways, as discussed in the section on finite size effects. 

This would imply that for sufficiently large system sizes one would find a more significant variation in the energies 
across the Brillouin zone, as only some of the states could be obtained from the tower. One would also find for 36 site 
systems that the spinon solenoid fluxes would change under a twist in spin boundary conditions, and so the fermionic 
states at different k would lose their degeneracy, leading to a change in energy of some of the k ^ states when 
varying boundary conditions. 

Another possibility is that, even for infinite systems, short-distance effects lead to the production of a staggered- 
mass pattern, in the style of the "perfect hexagon" state |H| discussed above. The staggered mass pattern gives rises 
to an enlarged unit cell, and makes it possible to get low energy states at the other points in the Brillouin zone. 



X. FURTHER COMPARISON WITH NUMERICS 



In addition to the existence of the low energy states, we make further comparison with numerics for the many-body 
density of states and the dimer-dimer correlation functions. 

Assuming the existence of a low-energy bosonic mode with linear dispersion relation, so that the single-particle 
density of states scales linearly with energy, one would expect the many-body density of states at energy E to scale 
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for large systems as an exponential of E 2 / 3 . The quadratic behavior of experimentally measured specific heat is in 
agreement with this. 

Since we have argued that the low-energy states in numerical calculations are largely "tower" states, it is impossible 
to extract the density of states in a large system from the finite size density of states. The true exponential growth 
can only be seen when the k ^ modes become important. 

In this regard, the numerically measured [pit quadratic many-body density of states at very low energies does not 
say anything about the dispersion of the Goldstone modes. Instead, it is a reflection of the fact that if a finite number 
of "tower" modes are excited then the many-body density of states is a power law. 

For the dimer-dimer calculations we can make a more direct comparison with a numerical calculation of these 
correlations on a 36 site system Q. In the thermodynamic limit, the dimer-dimer correlation function should be 
long-ranged, reflecting the existence of a non-zero TOi2- 

However, there is also a short-range fermionic contribution to the dimer-dimer correlation functions, and for 36 
site systems, so that the system size is smaller than the correlation area of the fermions, we can ignore the effect 
of a non-zero ra\2 on the dimer-dimer correlation function and directly study the correlation functions of massless 
fermions. 

If we are interested in a dimer-dimer correlation function 

C (<ii)(M) - ((Si ■ S,)(S k ■ S t )) - ((Si ■ S d ))((S k ■ Si)) (46) 

we can express this, under the assumption of weak spinon interaction, directly in terms of the spinon Green's functions. 
Writing each spin operator in terms of spinons and considering various contractions we obtain 

C(i,j)(fc,0 = —^ e (GijGjkGkiGu) — 4Re(GijGjiGikGki) — Re(GikGkjGjiGu) + -^\Gik\ 2 \Gji\ 2 + -^\Gij\ 2 \Gji\ 2 (47) 

To obtain quantitatively accurate answers, we must include the effects of projection within an approximation like that 
used above, projecting on site i, j, fc, I, requiring that there be one fermion on each of these four sites. We have done 
this at 3 different levels of approximation. 

At the lowest level, we have noted that the wavefunction before projection will have one fermion on each of these 
sites is roughly the product of the probability that it will have one fermion on each of i, j by the probability that it 
will have one fermion on each of k,l. Given that | G^j | = 0.221383... for neighboring i, j, we should replace equation 
©by 

G(i,j)(kj) — K*(—4:Re(GijGjkGkiGii) — 4Re(GijGjiGikGki) — Re(GikGkjGjiGn) + -|Gifc| 2 |G.//| 2 + — IG^- 1 2 1 2 ^ 

(48) 

where 

K = (1.9264) 2 (49) 

is the desired factor of ( 1+1 g| G |4 ) 2 - 

At a more refined level, we have carried out the computation projecting on all four sites exactly. At the third level 
of approximation we have started to project out onto additional sites as well. We have calculated the correlation 
functions in this approximation, in the limit of an infinite system size, for pairs of bonds that can both be written 
in the same 12-site unit cell. We compare the result to results from numerics on a 36-site system Q. We do not 
consider pairs of bonds that cannot be written in the same unit cell, as at this separation, finite-size effects will 
become important in the numerics and the comparison will become impossible. We could in principle improve on our 
comparison with numerics by computing the spinon correlation functions in a 36-site system also, in which case it 
should be possible to compare all bonds, but we have not done this. 

The results are shown in table I, where in the column Theory I we show the simplest approximation, and in the 
column Theory II we show the second level of approximation. Qualitatively the theory works quite well on the signs 
even at this level, getting 9 out of 12 correct. The only signs that the theory gets wrong at this level occur when the 
theory predicts a very small value (< .01) for the correlation function. 

To improve this result, we included the third level of approximation in which we also project out onto a fifth site 
(to) for those dimer correlation functions such that there is one and only one site (to) which neighbors both dimers. 
In the column Theory III we show this level of approximation, as well as the particular site (to) that we picked. Once 
this is done all the signs work out for 11 out of 12 correlations, and the qualitative agreement is almost perfect. 

The magnitudes work out less well, as most of the dimer-dimer correlations are far too small within the spinon 
calculation. However, RVB calculations are quite poor at getting long range correlations without including some 
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gauge fluctuations. For example, in the one-dimensional Heisenberg antiferromagnet, the spins on the same sublattice 
are uncorrelated j37|. By including gauge fluctuations, this result can be substantially improved |]3^ , ^9| . 

For the largest positive correlation function, CVg 7)(i,8)j the magnitude does work out well even at the simplest 
approximation. However, for the largest negative function, Cy6,7)(h 5)) the magnitude is off by roughly a factor of 
5, until we go to the third level of approximation, at which point the magnitude becomes roughly correct. We can 
hope that a better inclusion of fluctuations will improve these results, just as it has done for the one-dimensional 
chain. Perhaps projecting on an entire 12-sitc cell would give better results, as suggested by the improvement in the 
results in column III. Detailed calculations of projection for some trial wavefunctions on the kagome lattice have been 
performed by Hsu and Schofield [p~5[. A similarly detailed calculation for our parent state would be of interest. 



XI. CONCLUSION 



In conclusion, we have constructed an RVB state on the kagome lattice, the "parent state", which has a Dirac 
structure. Consideration of the various mass perturbations to the Dirac equation unifies several other previously 
suggested long-range states of the kagome lattice Heisenberg antiferromagnet, with the exception of the BCS |ll| 
state and the bosonic Sp(2N) state [ p~4j| . 

While at the projected mean-field level the chiral spin liquid appears to be the best RVB state, we have argued 
by a renormalization group treatment that fluctuations provide a mechanism for stabilizing a state with a nonchiral 
mass term. The numerical evidence also argues against a chiral state. We have then proceeded to explicit comparison 
with numerics, taking as input only one quantity, the triplet mass gap. 

The physical idea behind our construction is that, given massless fermions, the system must ultimately try to break 
some symmetry to give mass to the fermions. There are many ways of doing this, but they all correspond to either 
introducing chirality, or to introducing some kind of spin-solid state in which the system dimerizes the tij . If we ignore 
the chiral state, we must have some kind of spin solid: we have proposed one possible spin solid and its attendant 
pseudo-Goldstone excitations. Other spin solids are possible, and may in fact be realized in favor of our proposal, but 
the general principle that the massless fermions must break a symmetry and produce pseudo-Goldstone excitations 
should be robust. Two other possibilities for spin solids that must be considered are the constant M 6 solid, and the 
perfect hexagon system with spatially varying Mq. 

Experimentally, it may be possible in principle to detect the spin-Peierls order. Although this as appears as long- 
range order only in 4 spin correlation functions, it should give rise to a short-range oscillatory piece in the spin-spin 
correlation function. Since the spin-spin correlation function decays exponentially and the dimerization is weak, this 
would be very difficult to detect, but in principle it is possible via neutron scattering. 

Theoretically, more work is needed on the fluctuations about the parent state. Doing a Gutzwiller projection of 
the wavefunctions on a 36-site lattice should enable much more direct comparison with numerics, especially for the 
dimer-dimer correlation functions. 

Numerically, it may be possible to confirm the identification of the low-energy states with tower states from sym- 
metry breaking. If one computes a dimer-dimer correlation function with the same pair of sites taken at two different 
times it may be possible to see the oscillations in the mass field. This may be difficult, though, given the relatively 
weak amount of dimerization present and the problem of extracting the contribution of the mass term to dimer-dimer 
correlation functions from the background of the fermionic contribution. Similarly, one can try to compute a suscep- 
tibility to spin-Peierls ordering by explicitly dimerizing coupling constants J in equation (Q). If indeed the system 
wishes to spontaneously order in the thermodynamic limit, then the susceptibility to dimerization should be large. 
This may make it possible to unambiguously determine whether the system prefers the M12 or other ordering pattern. 

It might also be interesting numerically to look at model systems in which the Hamiltonian has additional terms 
coupling to the chirality operator on each triangle. This might make it possible to probe the stability of the system 
to the chiral mass term, as well as providing some interesting states in which the Dirac fermions are moving in a large 
net magnetic field. 



XII. APPENDIX: AN INTERESTING FLAT BAND CASE 



When the system has flux 7r/4 through each triangle and flux it/ 2 through each hexagon the band structure becomes 
very peculiar. Using a 12-site unit cell, we find that the lowest band is doubly degenerate, and almost exactly flat. 
The next band above that is quadruply degenerate and exactly flat. The higher bands, which are all empty, are not 
flat. 
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It is very unlikely that any such state could be stabilized. We have argued above the the kagome lattice is not a 
chiral spin liquid. However, it may be possible to add chirality operators to the Hamiltonian to tune the flux through 
the triangles. 

In this case, the physics of the flat band state would be very amusing. We are used to the fact that in, for example, 
the one-dimensional Heisenberg antiferromagnet, one can deduce that the spin-1 excitations are composite objects of 
two spinons by looking at the excitation spectrum. Since the energy of the spin-1 object is a sum of two different 
energies, there is a continuum of possible energies. When, however, one of the two spinons is a hole excitation from a 
flat band, the energy is constant over the band, and there is no sign of the composite nature of the spin waves when 
looking at their energy spectrum. One would have a situation with one spinon hopping freely over the lattice, while 
the other spinon sits unmoving on a given unit cell! 
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TABLE I. Comparison of C^^kj) between numerical and spinon calculations. See figure 1 for labeling of points in unit 
cell. See text for discussion of various approximations. 



(i,j)(k,l) 


Numerics 


Theory I 


Theory II 


Theory III 


(m) 


(6,7)(1,8) 


0.04337 


0.08242 


0.06249 






(6,7)(8,2) 


-0.01416 


-0.01107 


-0.001318 


-0.01686 


(1) 


(6,7)(9,3) 


-0.00646 


-0.001041 


-0.0004279 


-0.022153 


(8) 


(6,7)(10,3) 


0.01178 


0.001043 


0.0005004 






(6,7)(10,4) 


-0.01045 


-0.001040 


-0.0005008 






(6,7)(11,5) 


-0.06510 


-0.02499 


-0.0112132 


-0.083829 


(12) 


(11,10)(7,8) 


0.01221 


0.018885 


0.00938 






(11,10)(9,8) 


-0.00113 


-0.01824 


-0.0025 


-0.01395 


(3) 


(11,10)(2,9) 


0.00108 


0.0001619 


0.00559 


-0.007171 


(3) 


(6,7)(11,4) 


0.01322 


-0.0009778 


-0.000327 


0.001416 


(12) 


(11,10)(2,8) 


0.00045 


-0.004377 


-0.001964 


0.00386 


(9) 


(6,7)(9,2) 


0.01322 


-0.0009788 


-0.000328 


0.001415 


(8) 
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4 

FIG. 1. 12-site unit cell, with fluxes indicated. Small numbers are used to label points for reference. 




FIG. 2. Band structure for the parent state. We scan along varying momenta in the x direction, at vanishing momentum in 
the f + direction. 
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FIG. 3. Mass perturbation to produce Mi . Bold lines are increased in strength, dotted lines are reduced in strength. 






FIG. 5. 7r-flux phase on square lattice with odd number of sites. Numbers label different squares, each containing 7r ± f flux. 





+ 


+ 





FIG. 6. Enlarged torus to compute wavefunctions. 
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